Randomizing world trade. I. A binary network analysis 
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The international trade network (ITN) has received renewed multidisciplinary interest due to 
recent advances in network theory. However, it is still unclear whether a network approach conveys 
additional, nontrivial information with respect to traditional international-economics analyses that 
describe world trade only in terms of local (first-order) properties. In this and in a companion paper, 
we employ a recently proposed randomization method to assess in detail the role that local properties 
have in shaping higher-order patterns of the ITN in all its possible representations (binary/weighted, 
directed/undirected, aggregated/disaggregated by commodity) and across several years. Here we 
show that, remarkably, the properties of all binary projections of the network can be completely 
traced back to the degree sequence, which is therefore maximally informative. Our results imply 
that explaining the observed degree sequence of the ITN, which has not received particular attention 
in economic theory, should instead become one the main focuses of models of trade. 



PACS numbers: 89.65. Gh; 89.70.Cf; 89.75.-k; 02.70.Rr 

I. INTRODUCTION 

The network of import/export trade relationships 
among all world countries, known in the literature as the 
International Trade Network (ITN) or the World Trade 
Web (WTW), has received a renewed multidisciplinary 
interest in recent years jD-ll3| , due to impressive advances 
in both empirical and theoretical approaches to the study 
of complex networks [T4lfT6] . A number of robust pat- 
terns in the structure of this network have been empir- 
ically observed, both in its binary (when only the pres- 
ence of a trade interaction is considered, irrespective of 
its intensity) and weighted (when also the magnitude of 
trade flows is taken into account) description. These styl- 
ized facts include local properties as well as higher-order 
patterns. Local properties involve direct (first-order) in- 
teractions alone, resulting in simple quantities such as 
node degree (the number of trade partners of a country) , 
node strength (total trade volume of a country) , and their 
directed-network analogues (i.e., when these statistics are 
computed taking into account edge/trade directionality). 
Higher-order characteristics are more complicated struc- 
tural properties that also involve indirect interactions, i.e. 
topological paths connecting a country to the neighbors 
of its neighbors, or to countries farther apart. Exam- 
ples include degree-degree correlations, average nearest- 
neighbor indicators, and clustering coefficients, to name 
just a few of them. 

In general, local and higher-order topological proper- 
ties are not independent of each other. In particular, 



even if one assumes that the network is formed as the 
result of local constraints alone, with higher-order prop- 
erties being only the mere outcome of a specification of 
these constraints, it turns out that so-called structural 
correlations are automatically generated. Structural cor- 
relations sometimes appear as complicated patterns that 
might be confused with genuine correlations involving 
higher-order statistics, and interpreted as the presence 
of an additional level of topological organization. There- 
fore, in any real network it is important to characterize 
structural correlations and filter them out in order to as- 
sess whether nontrivial effects due to indirect interactions 
are indeed present. 

In the specific case of the ITN, this problem is partic- 
ularly important to assess whether the network formal- 
ism is really conveying additional, nontrivial information 
with respect to traditional international-economics anal- 
yses, which instead explain the empirical properties of 
trade in terms of country-specific macroeconomic vari- 
ables alone. Indeed, the standard economic approach 
to the empirics of international trade |17j has tradition- 
ally focused its analyses on the statistical properties of 
country-specific indicators like total trade, trade open- 
ness (ratio of total trade to GDP, i.e. Gross Domestic 
Product), number of trade partners, etc., that can be 
easily mapped to what, in the jargon of network anal- 
ysis, one denotes as local properties or first-order node 
characteristics. Ultimately, understanding whether net- 
work analyses go a step beyond with respect to stan- 
dard trade theory amounts to assess the effects of in- 
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direct interactions in the world trade system. Indeed, a 
wealth of results about the analysis of international trade 
have already been derived in the macroeconomics liter- 
ature [17) without making explicit use of the network 
description, and focusing on the above country-specific 
quantities alone. Whether more recent analyses of trade, 
directly inspired by the network paradigm jD-112]. are 
indeed conveying additional and nontrivial information 
about the structure of international import/export flows, 
crucially depends on the answer to the above question. 
Some network-inspired studies have already tried to ad- 
dress this problem, but with ambiguous results. In some 
cases, it was suggested that local properties are enough 
to explain higher-order patterns (3] [18] , while in others 
the opposite conclusion was reached [TO] . However, pre- 
vious analyses of the ITN focused on heterogeneous rep- 
resentations (either binary [2] [3] or weighted [6] [TOJ ES] , 
either directed [U [20] or undirected [2] [3] , either ag- 
gregated [H [3] [TO] or disaggregated [TO] in separate com- 
modities) and using different datasets, making consistent 
conclusions impossible. 

In this and in a companion paper |21| . we explicitly 
address this problem and exploit a recently proposed an- 
alytical method [55] to obtain, for any given topological 
property of interest, the value of the corresponding quan- 
tity averaged over the family of all randomized variants of 
the ITN that preserve the observed local properties. This 
allows us to identify empirical deviations from locally- 
induced structural correlations. Null models are used in 
our exercises to uncover significant features of the net- 
work and to understand to which extent some network 
statistics are sufficient to explain other network statis- 
tics. Our analysis is not however involved in explaining 
the underlying causal mechanisms shaping the network. 
Therefore, throughout this and its companion paper, we 
shall use the term "explaining" in a weak sense. For 
example, finding that a local network statistics X "ex- 
plains" a higher-order network statistics Y in our null 
model will signal the presence of a strong correlation be- 
tween the two statistics, so that X can be sufficient to 
fully reproduce Y in the network. Of course, we do not 
aim at using our null model to identify subtle causal links 
between X and Y, which in the real- world may be caused 
e.g. by some omitted variables that cause in a proper 
way the high observed correlation between X and Y. 

In this first paper, we focus on the ITN as a binary 
network. We find that higher-order patterns of all bi- 
nary (either directed or undirected) projections of the 
ITN are remarkably well explained by local properties 
alone (the degree sequences) . This result is robust to dif- 
ferent levels of commodity aggregation: even if with an 
increasing scatter, the degree sequence preserves its com- 
plete informativeness as more disaggregated and sparser 
commodity-specific networks are considered. Moreover, 
we perform a temporal analysis and check the robust- 
ness of these results over time. Therefore we obtain, for 
the first time in this type of study, a detailed and homo- 
geneous assessment of the role of local properties across 



different representations of the trade network, using var- 
ious levels of commodity aggregation, and over several 
years. From an international-trade perspective, our re- 
sults indicate that binary network descriptions of trade 
can be significantly simplified by considering the degree 
scqucncc(s) only. In other words, in any binary repre- 
sentation of the ITN, the degree sequence turns out to 
be maximally informative, since its knowledge conveys 
almost the entire information about the topology of the 
network. 

In the companion paper [21] . we show that the pic- 
ture changes completely when considering the ITN as a 
weighted network. We find that the ITN is an excel- 
lent example of a network whose local topological prop- 
erties cannot be deduced from its local weighted prop- 
erties. These results highlight an important limitation 
of current economic models of trade, that do not aim at 
explaining or reproducing the observed degree sequence 
but focus more on the structure of weights |23j . In other 
words, standard models of trade in economics have been 
focusing only on explaining the (positive) flow between 
any two countries, disregarding to a great extent theories 
that are able to account for the determinants of the cre- 
ation of a link (i.e. the transition from a zero trade-flow 
to a positive trade flow) . The observed extreme informa- 
tiveness of the degree sequence leads us to conclude that 
such models should be substantially revised in order to 
explicitly include the degree sequence of the ITN among 
the key properties to reproduce. 

II. DATA AND METHODS 

This Section describes the data we use to construct the 
various representations of the network in this and in the 
following paper [5T], discusses how the country-specific 
properties that are usually considered in world-trade eco- 
nomics translate into local topological properties of the 
ITN, and discusses how these properties should be kept 
as constraints of our analysis using an appropriate net- 
work randomization method. 



A. The International Trade Network 

We use yearly bilateral data on exports and imports 
from the United Nations Commodity Trade Database 
(UN COMTRADE) 04] from year 1992 to 2002. We 
have chosen this database because, despite its relatively 
short time interval (11 years), it contains trade data be- 
tween countries disaggregated across commodity cate- 
gories. This allows us to perform our analyses both at the 
aggregate level (total trade flows) and at the commodity- 
specific level, e.g. investigating whether local properties 
are sufficient to explain higher-order ones in commodity- 
specific networks of trade. 

In order to perform a temporal analysis and allow com- 
parisons across different years, we restrict ourselves to a 
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balanced panel of N = 162 countries that are present in 
the data throughout the time interval considered. As to 
the level of disaggregation, we choose the classification 
of trade values into C = 97 possible commodities listed 
according to the Harmonized System 1996 (HS1996) [15] . 
Accordingly, for a given year t we consider the trade value 
e?j (t) corresponding to exports of the particular commod- 
ity c (c = 1, . . . C). Since, for every commodity, exports 
from country i to country j are reported twice (by both 
the importer and the exporter) and the two figures do not 
always match, we follow Ref. |13] and only employ the 
flow as reported by the importer. Besides commodity- 
specific data, we also compute the total value e^(t) of 
exports from country i to country j as the sum over the 
exports of all C = 97 commodity classes: 



c=l 



(1) 



The particular aggregation procedure described above, 
which coincides with the one performed in Ref. |13j . al- 
lows us to compare our analysis of the C commodity- 
specific networks with a (C + l)-th aggregate network, 
avoiding possible inconsistencies between aggregated and 
disaggregated trade data. We stress that the resulting ag- 
gregated network data are in general different from those 
used in other analyses [31 El H] of the same network. 
Nonetheless, as we show below, when we analyze net- 
work properties that have also been studied in previous 
studies of aggregate trade, we find perfect agreement. 

The quantities {e£,-(i)} (where c = 0, . . . C) defined 
above are the fundamental data that allow us to obtain 
different possible representations of the trade network, as 
well as the corresponding randomized counterparts (see 
below for the units of measure we adopted). When we 
regard the ITN as a weighted directed network, we define 
the weight of the link from country i to country j in year 
t for commodity c as 



c = 0, 



(2) 



where <E N denotes the nearest integer to the non- 
negative real number x. When we adopt a weighted 
but undirected (symmetrized) description, we define the 
weight of the link between countries i and j in year t for 
commodity c as 



/(*) 
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(t) 



e?,-(*)+e£(t) 



c = 0, 



C (3) 



Therefore, in both the directed and undirected case, 
wfj (t) is an integer quantity. Since in both cases we shall 
be interested in tracking the temporal evolution of most 
quantities, we also define rescaled weights (relative to the 
total yearly trade flow) as follows: 



(*) 



t(t) 



0,...C 



(4) 



where in the directed case w^{t) is given by Eq. (|2j) 
and Wt ot {t) = J2iY^j J ti w ij('t'j (the double sum runs 
over all N(N — 1) ordered pairs of vertices), while 
in the undirected case wfj(t) is given by Eq. (|3j) and 
w% ot (t) = J2i ^j<i (t ne double sum runs over 

all the N(N — l)/2 unordered pairs). In such a way, 
trend effects are washed away and we obtain adimen- 
sional weights that are automatically deflated, allowing 
consistent comparisons across different years and differ- 
ent commodities. 

In the binary representations of the network, we draw a 
link from i to j whenever the corresponding weight wfj is 
strictly positive. If Q(x) denotes the step function (equal 
to 1 if x > and otherwise), the adjacency matrix of the 
binary projection of the network in year t for commodity 
c is 



aut) = e[wut)] c = 0,...C 



(5) 



where tOy(i) is given either by Eq. (|2j) or by Eq. (|3j), 
depending on whether one is interestea in a directed or 
undirected binary projection of the network respectively. 

For each of the C + 1 commodity categories, we 
can consider four network representations (binary undi- 
rected, binary directed, weighted undirected, weighted 
directed). When reporting our results, we will first 
describe the aggregated networks (c = 0) and then 
the disaggregated (commodity-specific) ones. In par- 
ticular, among the 97 commodity classes, we will fo- 
cus on the 14 particularly relevant commodities identi- 
fied in Ref. |13j . which are reported in table |TJ These 
14 commodities include the 10 most traded commodi- 
ties (c = 84,85,27,87,90,39,29,30,72,71 according to 
the HS1996) in terms of total trade value (following the 
ranking in year 2003 [32]), plus 4 classes (c = 10, 52, 9, 93 
according to the HS1996) which are less traded but more 
relevant in economic terms. Taken together, the 10 most 
traded commodities account for 56% of total world trade 
in 2003; moreover, they also feature the largest values 
of trade value per link (also shown in the table). The 
14 commodities considered account together for 57% of 
world trade in 2003. As an intermediate level of aggrega- 
tion, we shall also consider the networks formed by the 
sum of these 14 commodities. The original data {e^-(t)} 
are available in current U.S. dollars (USD) for all com- 
modities; however, due to the different trade volumes 
involved, we use different units of measure for different 
levels of aggregation [4"B] . 



B. Controlling for local properties 

As we mentioned, our main interest in the present work 
is assessing whether higher-order properties of the ITN 
can be simply traced back to local properties, which are 
the main focus of traditional macroeconomic analyses of 
international trade. Such standard country-specific prop- 
erties include: total exports, total imports, total trade 
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HS Code Commodity 


Value (USD) Value per link (USD) % of aggregate trade 


84 


Nuclear reactors, boilers, machinery and 
mechanical appliances; parts thereof 


5.67 x 10 11 


6.17 x 


10' 


11.37% 


85 


Electric machinery, equipment and parts; 
sound equipment; television equipment 


5.58 x 10 11 


6.37 x 


10 7 


11.18% 


27 


Mineral fuels, mineral oils & products of 
their distillation; bitumin substances; min- 
eral wax 


4.45 x 10 11 


9.91 x 


10 7 


8.92% 


87 


Vehicles, (not railway, tramway, rolling 
stock); parts and accessories 


3.09 x 10 11 


4.76 x 


10 7 


6.19% 


90 


Optical, photographic, cinematographic, 
measuring, checking, precision, medical or 
surgical instruments/apparatus; parts & 
accessories 


1.78 x 10 11 


2.48 x 


10 7 


3.58% 


39 


Plastics and articles thereof. 


1.71 x 10 11 


2.33 x 


10 7 


3.44% 


29 


Organic chemicals 


1.67 x 10 11 


3.29 x 


10 7 


3.35% 


30 


Pharmaceutical products 


1.4 x 10 11 


2.59 x 


10 7 


2.81% 


72 


Iron and steel 


1.35 x 10 11 


2.77 x 


10 7 


2.70% 


71 


Pearls, precious stones, metals, coins, etc 


1.01 x 10 11 


2.41 x 


10 7 


2.02% 


10 


Cereals 


3.63 x 10 10 


1.28 x 


10 7 


0.73% 


52 


Cotton, including yarn and woven fabric 
thereof 


3.29 x 10 10 


6.96 x 


10 6 


0.66% 


9 


Coffee, tea, mate & spices 


1.28 x 10 10 


2.56 x 


10 6 


0.26% 


93 


Arms and ammunition, parts and acces- 
sories thereof 


4.31 x 10 9 


2.46 x 


10 6 


0.09% 


ALL 


Aggregate (all 97 commodities) 


4.99 x 10 12 


3.54 x 


10 s 


100.00% 



TABLE I: The 14 most relevant commodity classes (plus aggregate trade) in year 2003 and the corresponding total trade value 
(USD), trade value per link (USD), and share of world aggregate trade. From Ref. |13j . 

(sum of total exports and total imports) , trade openness 
(ratio of total trade to GDP), the number of countries 
whom a country exports to and imports from, the total 
number of trade partners (irrespective of whether they 
are importers or exporters, or both). All these quan- 
tities can be simply obtained as local sums over direct 
interactions (countries one step apart) in a suitable rep- 
resentation of the network. 

For instance, the number of trade partners of country i 
is simply the number of neighbors of node i in the binary 
undirected projection, i.e. the degree 

h = Y^ a l3 (6) 

In the above equation and in what follows, we drop the 
dependence of topological quantities on the particular 
year t for simplicity. We also drop the superscript c speci- 
fying a particular commodity, as all the formulas hold for 
any c. This means that, if the aggregated network of total 
trade is considered, then aij and Wij represent the aggre- 
gate quantities a®j and u>^ , where the commodity c = 
formally represents the sum over all commodities, as in 
Eq. ([I]). Otherwise, if the commodity-specific network 
involving only the trade of the particular commodity c 
(with c > 0) is considered, then and w,j represent the 
values a?j and w\a for that commodity. 

The number of countries whom a country exports to 
and imports from are simply the two directed analogues 



(the out-degree k° ut and the in-degree k™ 1 respectively) 
of the above quantity in the binary directed description: 

fcr* = E«Hi (?) 

kT = (8) 

Similarly, as evident from Eq. ([3]), country i's total 
trade coincides with twice the sum of weights reaching 
node i in the weighted undirected representation, i.e. the 
strength 

s i=^2 w ij ( 9 ) 

Finally, total exports (imports) of country i are simply 
the sum of out-going (in-coming) weights in the weighted 
directed representation of the ITN. These quantities are 
known as the out-strength s° ut and in-strength s™ of node 
i: 

s r = j>« (io) 

S T = 5>,i (ii) 

Another country-specific property which is widely used 
as an explanatory variable of trade patterns is the GDP 
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or the per capita GDP (i.e. the ratio of GDP to popula- 
tion). This property is sometimes used to rescale trade 
values, as in the case of trade openness which is defined 
as a country's ratio of total trade to GDP. Unlike the 
quantities discussed above, the GDP is not a topological 
entity. Nonetheless, it is empirically observed to be pos- 
itively (and strongly) correlated with the degree [5] and 
with node strength [12] (we will comment more on this in 
Section 111). Therefore, even if this is not the main aim 



of the present work, one should be aware that assessing 
the role of local topological properties also indirectly im- 
plies, to a large extent, assessing the role of the GDP of 
countries. 



denote (X) as a randomized property, and its value as the 
randomized value of X, even if technically no randomiza- 
tion process has been required (all the results have been 
obtained analytically). Similarly, we imagine the graph 
ensemble as a rewired version of the original network, 
even if no rewiring has taken place explicitly. 



III. THE ITN AS A BINARY UNDIRECTED 
NETWORK 



C. Rewiring the ITN 

We showed that, in a network language, the standard 
country-specific properties used to characterize world 
trade translate into simple local topological properties of 
the ITN. This naturally implies that, in our analysis, it is 
important to consider a null model of the ITN where such 
properties are enforced as constraints, and the topology 
is otherwise maximally random. Different methods that 
produce randomized ensembles of networks with given 
constraints exist [22j [24V 33J . As we mentioned, we aim 
at studying many topological properties of several differ- 
ent representations and temporal snapshots of the ITN. 
Therefore, we need a fast method that can deal with 
many networks in a relatively short time, and treat bi- 
nary, weighted, directed and undirected graphs in a con- 
sistent fashion. To this end, we employ the maximum- 
likelihood method introduced in ref . [22] , which provides 
the expectation values (over the randomized ensemble) of 
the desired topological properties analytically, in contrast 
with alternative methods [31, 32 which require to explic- 
itly generate many randomized variants of the real net- 
work computationally. Moreover, the method is density- 
independent and works for both sparse and dense net- 
works. By contrast, other (analytical or computational) 
approaches are density-dependent and not optimized for 
dense networks: the Chung-Lu (analytical) approach [33] 
works only for sparse networks, and the Maslov-Sneppen 
(computational) algorithm [311122] becomes too time con- 
suming for dense networks. Since the ITN is an unusu- 
ally dense network, the maximum-likelihood method is 
the natural choice that allows us to perform a detailed 
analysis, covering all possible representations across sev- 
eral years, which would otherwise require an impressive 
amount of time. 

In the Appendix we describe the maximum-likelihood 
method in some detail, in particular its application to 
the topological properties of interest for the present case 
study. Given any topological property X, the method 
provides the average value (X) of X across the ensem- 
ble of random graphs with the same average (across the 
ensemble itself) constraints as the real network. For sim- 
plicity, in this and in the companion paper we sometimes 



As we mentioned in Section |II A[ in its binary rep- 
resentation the ITN is defined as a graph whose edges 
report the presence of trade relationships among world 
countries, irrespective of the intensity of these relation- 
ships. The binary representation of the ITN can be either 
undirected or directed, depending on whether one is in- 
terested in specifying the orientation of trade flows. In 
both cases, the complete information about the topology 
of the network is encoded in the adjacency matrix A, 
whose entries {fly} are defined as in Eq. Q. 

In the simplest case, the presence of at least one of 
the two possible trade relationships between any two 
countries i and j (either from i to j or from j to i) 
is represented as one undirected edge between nodes i 



and j. Therefore 



and A is a symmetric ma- 



trix. In this binary undirected description, as shown in 
Eq. Q, the local constraints {C a } are the degrees of 
all vertices, i.e. the degree sequence {ki}. Therefore, the 
maximum-likelihood randomization method |22j (see Ap- 
pendix) works by specifying the constraints {C a } = {ki} 
and allows us to write down the probability of any graph 
G in the grandcanonical ensemble, which is uniquely 
specified by its generic adjacency matrix A. As sum- 
marized in Appendix [B] this allows us to easily obtain 
the expectation value (X) , formally defined in Eq. ( A7 1 , 
of any property X across the ensemble of binary undi- 
rected graphs whose expected degree sequence is equal to 
the empirical one. Note that, among the possible prop- 
erties, the degree of vertices plays a special role, as its 
expectation value (ki) is exactly equal to the empirical 
value ki, as required by the method. Therefore the values 
{ki} are useful control parameters and can be efficiently 
used as independent variables in terms of which other 
properties X can be visualized. 



For the sake of simplicity, in Sections HI A| and |III B 



we 

first report the results of this analysis on a single snap- 
shot of the commodity-aggregated network (the last year 
in our temporal window, i.e. 2002). Then, we discuss 
the robustness of our results through time by tracking 
them backwards in Section |HI C| We finally consider the 
disaggregated analysis of commodity-specific networks in 
Section [TlTDl 
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FIG. 1: (Color online) Average nearest neighbor degree k" n 
versus degree fc; in the 2002 snapshot of the real binary undi- 
rected ITN (red points), and corresponding average over the 
maximum-entropy ensemble with specified degrees (blue solid 
curve) . 
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FIG. 2: (Color online) Clustering coefficient C; versus degree 
ki in the 2002 snapshot of the real binary undirected ITN 
(red points), and corresponding average over the maximum- 
entropy ensemble with specified degrees (blue solid curve). 



A. Average nearest neighbor degree 

We start with the analysis of the aggregated version of 
the ITN, representing the trade of all commodities (c = 
in our notation). In the following formulas, the matrix 
A therefore denotes the aggregate matrix A , where we 
drop the superscript for brevity. As a first quantity, we 
consider the average nearest neighbor degree (ANND) of 
vertex i, defined as 



E 



ijj kj 



J2j^i Efc^j "'/".;'■ 



E 



(12) 



and measuring the average number of partners of the 
neighbors of a given node i. The above quantity in- 
volves indirect interactions of length two, as evidenced 
from the presence of terms of the type aijajk in the def- 
inition. Whether these 2-paths are a simple outcome of 
the concatenation of two independent edges can be in- 
spected by considering the correlation structure of the 
network, and in particular by plotting fc" n versus fc, . The 
result is shown in Fig. [T] We observe a decreasing trend, 
confirming what already found in previous studies em- 
ploying different datasets [21 IU [12]. This means that 
countries trading with highly connected countries have a 
few trade partners, whereas countries trading with poorly 
connected countries have many trade partners. This cor- 
relation profile, known as disassortativity, might signal 
an interesting pattern in the trade network. However, if 
we compare this trend with the one followed by the corre- 
sponding randomized quantity (fc™ n ) (see Appendix[B]for 
its expression), we find that the two behaviors coincide. 
This is an important effect of structural constraints in a 
dense network [34 : contrary to what naively expected 
[35], even in a network where links are drawn randomly 
between vertices with given heterogeneous degrees, the 
ANND is not constant. This means that the degree se- 
quence constrains the correlation structure, and that it is 



impossible to have a flat profile (kf n independent of ki) 
unless one forces the system to display it by introducing 
additional mechanisms (hence additional correlations of 
opposite sign). 



B. Clustering coefficient 

A similar result is found for the behavior of the clus- 
tering coefficient Cj, representing the fraction of pairs of 
neighbors of vertex i which are also neighbors of each 
other: 



Ej^i Efc^iJ a ij a jkdki 

Hh - 1) 

Ej^j Efc^ij a ij a jkaki 
Ej / - E/-v ... a ij a ik 



(13) 



The clustering coefficient is a measure of the fraction of 
potential triangles attached to i that are actually real- 
ized. This means that indirect interactions of length 
three, corresponding to products of the type aijajkaki 
entering Eq. (13), now come into play. Again, we find a 
decreasing trend of Cj as a function of ki (see Fig.[2|. This 
means that trade partners of highly connected countries 
are poorly interconnected, whereas partners of poorly 
connected countries are highly interconnected. However, 
if this trend is compared with the one displayed by the 
randomized quantity (cj) (see AppendixjB]), we again find 
a very close agreement. This signals that in the ITN also 
the profile of the clustering coefficient is completely ex- 
plained by the constraint on the degree sequence, and 
does not imply the presence of meaningful indirect inter- 
actions on top of a concatenation of direct interactions 
alone. 

The above results show that the patterns observed in 
the binary undirected description of the ITN do not re- 
quire, besides the fact that different countries have spe- 
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cific numbers of trade partners, the presence of higher- 
order mechanisms as an additional explanation. On the 
other hand, the fact that the degrees alone are enough to 
explain higher-order network properties means that the 
degree sequence is an important structural pattern in its 
own. This highlights the importance of reproducing the 
observed degree sequence in models of trade. We will 
comment more about this point later on. 



C. Evolution of binary undirected properties 

We now check the robustness of the previous results 
through time. This amounts to perform the same anal- 
ysis on each of the 11 years in our time window rang- 
ing from 1992 to 2002. For each of these snapshots, we 
specify the degree sequence and evaluate the maximally 
random ensemble of binary undirected graphs. We then 
compare each observed property X with the correspond- 
ing average (X) (repeating the procedure described in 
Appendix [B]) over the null model for that specific year. 
We systematically find the same results described above 
for each and every snapshot. For visual purposes, rather 
than replicating the same plots shown above for all the 
years considered, we choose a more compact description 
of the observed patterns and portray its temporal evolu- 
tion in a simple way. As we now show, this also provides 
us with a characterization of various temporal trends dis- 
played by each topological property, conveying more in- 
formation than a fixed-year description of the trade sys- 
tem. 

We first consider the average nearest neighbor degree. 
For a given year, we focus on the two lists of vertex- 
specific values {kf n } and {(&"")} for the real and ran- 
domized network respectively. We compute the average 
(rrifcnn and mn.»»\) and the associated 95% confidence in- 
terval of both lists and plot them together as in Fig. [3^. 
We repeat this for all years and obtain a plot which in- 
forms us about the temporal evolution of the ANND in 
the real and randomized network separately. We find 
that the average value of the empirical ANND has been 
increasing steadily during the time period considered. 
However, the same is true for its randomized value, which 
is always consistent with the real one within the confi- 
dence intervals. This means that the null model com- 
pletely reproduces the temporal trend of degree-degree 
correlations. 

In principle, the increase of the ANND could be sim- 
ply due to an overall increase in link density. To fur- 
ther study this possibility, we have compared the yearly 
growth rate (X t /X t -i — 1) of the average ANND and of 
the link density in the period considered. We found two 
regimes: initially (from 1993 to 1997) the density has 
a larger (but decreasing) growth rate than the ANND, 
while from 1998 to 2002 onwards the two rates converge. 
Therefore it is useful to keep in mind that the evolution 
of the average ANND, as well as that of other average 
properties we consider below, is in general not merely 



reflecting the evolution of the overall link density. 

In Fig. [3)3 we also plot the temporal evolution of the 
standard deviations s^n and s/k nn ) (with associated 95% 
confidence intervals) of the two lists of values and 
{(fc™)}. We find that the variance of the empirical aver- 
age nearest neighbor degree has been decreasing in time, 
but once more this behavior is completely reproduced by 
the null model and therefore fully explained by the evo- 
lution of the degree sequence alone. Moreover, in Fig. [3b 
we show the Pearson (product-moment) correlation co- 
efficient rjnn j. (with 95% confidence interval) between 
{k™ n } and {ki}, and similarly the correlation coefficient 
r (k nn ),k between the randomized quantities {(&"")} and 
{ki} (recall that {(ki)} — {ki} by construction). This in- 
forms us in a compact way about the evolution of the de- 
pendence of the ANND on the degree, i.e. of the change 
in the structure of the scatter plot we showed previously 
in Fig. [T] We find that the disassortative character of 
the scatter plot results in a correlation coefficient close 
to — 1, which has remained remarkably stable in time 
across the interval considered, and always very close to 
the randomized value. 

The complete accordance between the real and ran- 
domized ANND in each and every snapshot is confirmed 
by Fig. [3fcl, where we show the correlation coefficient 
r fc"",( fc ""> ( w ith 95% confidence interval) between the em- 
pirical ANND, {fcf™}, and the randomized one, {(&"")}• 
We observe an approximately constant value close to 1, 
signaling perfect correlation between the two quantities. 
This exhaustively explains the accordance between the 
real and randomized ANND for all vertices, while the 
other three panels of Fig. [3] also inform about various 
overall temporal trends of the ANND, as we discussed. 

In Fig.|4]we show the same analysis for the values {ct} 
and {(cj)} of the clustering coefficient. In this case we 
observe an almost constant trend of the average cluster- 
ing coefficient (Fig. [4 1), a decreasing standard deviation 
(Fig. [4)3), and a stable strong anticorrelation between 
clustering and degree (Fig. [Ik) . Again, we find that the 
real and randomized values are always consistent with 
each other, so that the evolution of the empirical values 
is fully reproduced by the null model. This is confirmed 
by Fig. [4ji, which shows that the correlation between 
{ci} and {(ci)} is always very close to 1. As for the 
ANND, these results clearly indicate that the real and 
randomized values of the clustering coefficient of all ver- 
tices are always in perfect agreement, and that the tem- 
poral trends displayed by this quantity are completely 
explained by the evolution of the degree sequence. 



D. Commodity-specific binary undirected networks 

We complete our analysis of the ITN as a binary undi- 
rected network by studying whether the picture changes 
when one considers, rather than the network aggregat- 
ing the trade of all types of commodities, the individual 
networks formed by imports and exports of single com- 
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FIG. 3: (Color online) Temporal evolution of the properties 
of the nearest neighbor degree fc" n in the 1992-2002 snapshots 
of the real binary undirected ITN and of the corresponding 
maximum-entropy null model with specified degrees, a) av- 
erage of k" n across all vertices (red: real data; blue: null 
model, indistinguishable from real data), b) standard devi- 
ation of k" n across all vertices (red: real data; blue: null 
model, indistinguishable from real data), c) correlation co- 
efficient between fc" n and ki (red, upper symbols: real data; 
blue, lower symbols: null model), d) correlation coefficient 
between k" n and (k" n ). The 95% confidence intervals of all 
quantities are represented as vertical bars. 
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FIG. 4: (Color online) Temporal evolution of the proper- 
ties of the clustering coefficient d in the 1992-2002 snapshots 
of the real binary undirected ITN and of the corresponding 
maximum-entropy null model with specified degrees, a) aver- 
age of Ci across all vertices (red: real data; blue: null model, 
indistinguishable from real data), b) standard deviation of a 
across all vertices (red: real data; blue: null model, indistin- 
guishable from real data), c) correlation coefficient between 
d and ki (red, upper symbols: real data; blue, lower symbols: 
null model), d) correlation coefficient between Ci and (c;). 
The 95% confidence intervals of all quantities are represented 
as vertical bars. 



modities. To this end, we focus on the disaggregated 
data described in Section [II A| and we repeat the analysis 
reported above, by identifying the matrix A with various 
disaggregated matrices A c (with c > 0). 

We find that the results obtained in our aggregated 
study also hold for individual commodities. For brevity, 
we only report the scatter plots of the average near- 
est neighbor degree (Fig. [5| and clustering coefficient 
(Fig. [6]) for the 2002 snapshots of 6 commodity-specific 
networks. The 6 commodities are chosen among the top 
14 reported in Table [TJ In particular, we select the two 
least traded commodities in the set (c = 93, 9), two inter- 
mediate ones (c = 39, 90), the most traded one (c = 84), 
plus the network formed by combining all the top 14 
commodities, i.e. an intermediate level of aggregation 
between single commodities and the completely aggre- 
gated data (c = 0), which we already considered in the 
previous analysis (Figs. l]and[2|. With the addition of 
the latter, the results shown span 7 different cases or- 
dered by increasing trade intensity and level of commod- 
ity aggregation. Similar results hold also for the other 
commodities not shown. 

If we compare Fig. [5] with Fig. [I] we see that the 
trend displayed by ANND in the aggregated network is 
preserved, even if with a slightly increasing scatter, as 
sparser and less disaggregated commodity classes are con- 
sidered. Importantly, the accordance between real and 
randomized values is also preserved. The same is true 
for the clustering coefficient, cf. Fig. [6] and its compari- 
son with Fig. [2j These results indicate that the degree se- 



quence maintains its complete informativeness across dif- 
ferent levels of commodity resolution, and irrespective of 
the corresponding intensity of trade. Thus, remarkably, 
the knowledge of the number of trade partners involv- 
ing only a specific commodity still allows to reproduce 
the properties of the corresponding commodity-specific 
network. 

As a summary of our binary undirected analysis we 
conclude that, in order to explain the evolution of the 
ANND and clustering of the ITN, it is unnecessary to 
invoke additional mechanisms besides those accounting 
for the evolution of the degree sequence alone. Since 
the ANND and clustering already probe the effects of in- 
direct interactions of length two and three respectively, 
and since higher-order correlations involving longer topo- 
logical paths are built on these lower-level ones, the null 
model we considered here is very likely to fully reproduce 
the properties of the ITN at all orders. In other words, 
we found that in the binary undirected representation 
of the ITN the degree sequence is maximally informa- 
tive, as its knowledge allows to predict the higher-order 
topological properties of the network that we have ex- 
plored in this and in the companion paper. An inter- 
esting question is whether the degree sequence is also 
able to reproduce other higher-order network properties, 
such as path lengths, node centrality, etc.. Whereas this 
study does not explicitly address this question, we argue 
that the answer will be positive in the light of the very 
nature of the ITN. Its high density indeed implies that 
path lengths are almost never larger than 3. As a con- 
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FIG. 5: (Color online) Average nearest neighbor degree k" n 
versus degree fej in the 2002 snapshots of the commodity- 
specific (disaggregated) versions of the real binary undi- 
rected ITN (red points), and corresponding average over the 
maximum-entropy ensemble with specified degrees (blue solid 
curves), a) commodity 93; b) commodity 09; c) commodity 
39; d) commodity 90; e) commodity 84; f) aggregation of the 
top 14 commodities (see Table [I] for details). From a) to f), 
the intensity of trade and level of aggregation increases. 



FIG. 6: (Color online) Clustering coefficient Ci versus degree 
ki in the 2002 snapshots of the commodity-specific (disaggre- 
gated) versions of the real binary undirected ITN (red points) , 
and corresponding average over the maximum-entropy ensem- 
ble with specified degrees (blue solid curves), a) commodity 
93; b) commodity 09; c) commodity 39; d) commodity 90; 
e) commodity 84; f) aggregation of the top 14 commodities 
(see Table [I] for details). From a) to f), the intensity of trade 
and level of aggregation increases. 



sequence, network properties of order larger than 3 are 
typically well proxied by local properties. An example is 
the extremely high correlation between country centrality 
(measured e.g. in terms of betweenness centrality) and 
node degree, typically found in previous studies. The ro- 
bustness of this result across several years and different 
commodity classes strengthens our previous discussion 
about the importance of including the degree sequence 
among the focuses of theories and models of trade, which 
are instead currently oriented mainly at reproducing the 
weighted structure, rather than the topology of the ITN. 



IV. THE ITN AS A BINARY DIRECTED 
NETWORK 



We now consider the binary directed description of the 
ITN, with an interest in understanding whether the in- 
troduction of directionality changes the picture we have 
described so far. In the directed binary case, a graph G 
is completely specified by its adjacency matrix A which 
is in general not symmetric, and whose entries are a,j = 1 form a disaggregated analysis in Section IV D 



if a directed link from vertex i to vertex j is there, and 
a™ = otherwise. The local constraints {C a } are now 
the two sets of out-degrees and in-degrees of all vertices 
defined in Eqs.([7]) and ([8]), i.e. the out-degree sequence 
{k° ut } and the in- degree sequence {kf 1 }. In Appendixjc] 
we show how the randomization method enables in this 
case to obtain the expectation value (X) of a property X 
across the maximally random ensemble of binary directed 
graphs with in-degree and out-degree sequences equal to 
the observed ones. When inspecting the properties of 
the ITN and its randomized variants, the useful indepen- 
dent variables are now the values {k° ut } and {fc™} (or 
combinations of them) , since they are the special quanti- 
ties X whose expected value (X) coincides with the ob- 
served one by construction. Again, we first consider the 
2002 snapshot of the completely aggregated ITN (Sec- 
tions 



IV A and IV B I, then track the temporal evolution 



of the results backwards (Section IV C), and finally per- 
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A. Directed average nearest neighbor degrees 

We start with the analysis of the binary directed trade 
network aggregated over all commodities (c = 0). There- 
fore, in the following formulas, we set A = A . The 
average nearest neighbor degree of a vertex in a directed 
graph can be generalized in four ways from its undirected 
analogue. We thus obtain the quantities 
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(17) 



In the above expressions, indirect interactions due to 
the concatenation of pairs of edges are taken into ac- 
count according to their directionality, as clear from the 
presence of products of the type a^aki- A fifth possibil- 
ity is an aggregated measure based on the total degree 
k tot = k in + k out of vertices: 



k 



tot/tot 



k 



(18) 



The latter is a useful one to start with, as it provides a 
simpler analogue to the undirected quantity k™ we have 
already studied. At the same time, it must be noted that 
the two quantities are not trivially related since the total 
directed properties k\ ot and fc*°*/* ot carry more informa- 
tion than the corresponding undirected ones ki and fc" n , 
the difference being the local reciprocity structure of the 
network [20]. To see this, note that k\ ot = ki + k^, 
where k^ = a ij a ji (if a ij is t ne adjacency ma- 

trix of the directed network) is the reciprocated degree 
of vertex i, defined as the number of bidirectional links 
reaching i 20 , 36 , 37 . This quantity represents the num- 
ber of trade partners, acting simultaneously as importers 
and exporters, of country i. Therefore, studying total 
directed quantities also allows to assess whether the reci- 
procity structure of the directed network changes the 
results obtained in the undirected case (similar consid- 
erations apply to the directed clustering coefficients we 
introduce below). 

In Fig. [7] we plot \~-° l l tot as a function of k\ ot for the 
2002 snapshot of the binary directed ITN. The trend 
shown does not differ substantially from its undirected 
counterpart we observed in Fig. [T] In particular, we ob- 
tain a similar disassortative character of the correlation 
profile. Importantly, we find again a good agreement 
between the empirical quantity and its expected value 
(/c' ot ^* ot ) under the null model (obtained as in Appendix 



[C]).InFig.[| we show a more refined analysis by consider- 
ing all the four directed versions of the ANND defined in 
Eqs.([l"4"])-|l"7]), as well as their expected values under the 
null modeljsee Appendix[C]). We immediately see that all 
quantities still display a disassortative trend, with some 
differences in the ranges of observed values. Again, all 
the four empirical behaviors are in striking accordance 
with the null model, as the randomized curves (obtained 
as in Appendix [C]) show. This means that both the de- 
creasing trends and the ranges of values displayed by all 
quantities are well reproduced by a collection of random 
graphs with the same average in-degrees and out-degrees 
as the real network. 



B. Directed clustering coefficients 

We now consider the directed counterparts of the clus- 
tering coefficient defined in Eq. ( p~3| ) . Again, there are 
four possible generalizations depending on whether the 
directed triangles involved are of the inward, outward, 
cyclic or middleman type [38] : 
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The directed clustering coefficients are determined by in- 
direct interactions of length 3 according to their direc- 
tionality, appearing as products of the type aija k ia mn in 
the above formulas. At the same time, since they always 
focus on three vertices only, they capture the local occur- 
rence of particular network motifs [39] of order 3. A fifth 
aggregated measure, based on all possible directions, is 



k ^ij( a ij + a ji 

2[kl ot {k\ ot - 



)(a jk + a kj )(a ki + a ik ) 



2k 



(23) 

As for fc' ot ^*°*, the latter definition is a good starting 
point for a comparison with the undirected case. In Fig. [9] 
we show c* ot and (c- ot ) (see Appendix [c| as a function 
k\ ot for our usual snapshot. We see no fundamental dif- 
ference with respect to Fig. [2] Again, the randomized 
quantity does not deviate significantly from the empiri- 
cal one. 

We now turn to the four directed clustering coefficients 
defined in Eqs.( 19 1-(22 1. We show these quantities in 
10 as functions of different combinations of k™ and 
depending on the particular definition. As for the 



Fig. 

k° ut 



directed ANND, we observe some variability in the range 
of observed clustering values. However, all the quantities 
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FIG. 7: (Color online) Total average nearest neighbor de- 
gree /j* '/* * versus total degree k\ ot in the 2002 snapshot of 
the real binary directed ITN (red points), and corresponding 
average over the maximum-entropy ensemble with specified 
out-degrees and in-degrees (blue solid curve). 
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FIG. 9: (Color online) Total clustering coefficient c* 0< versus 
total degree k\ ot in the 2002 snapshot of the real binary di- 
rected ITN (red points), and corresponding average over the 
maximum-entropy ensemble with specified out-degrees and 
in-degrees (blue solid curve). 



are again in accordance with the expected ones under the 
null model (see Appendix [C]). 



C. Evolution of binary directed properties 

We now track the temporal evolution of the above re- 
sults by performing, for each year in our time window, 
an analysis similar to that reported in sec jlll C| for the 
undirected case. 

We start by showing the evolution of the total aver- 
age nearest neighbor degree m ^ ne f our panels of 



Fig. [TTJ where we plot the same properties considered 
previously for the undirected ANND in Fig. [3] We find 
that the temporal evolut ion of the average (Fig.[lljj,) and 



standard deviation (Fig. 11 d) of fc*°*/* ot j s essentially the 
same as that of the undirected fef™, apart from differences 
in the range of values. Similarly, th e co rrelation coeffi- 
cients between £*°'/* ot a nd /j* * (Fig. 11:), (fc*°*/ tot ) and 
(jfeft) = k\ ot (Fig. EL), k l ° t/tot and {kf t/tot ) (Fig. Illk) 



mimic their undirected counterparts, confirming that the 
perfect accordance between fc*°*/ tot anc ; (£;' ot / tot ) is stable 
over time, and that the disassortative trend of /j' */* * a s 
a function of k\ ot (Fig. [7]) is always completely explained 




20 40 60 80 100 120 20 40 60 80 100 120 

fcoot £001 

FIG. 8: (Color online) Directed average nearest neighbor de- 
grees versus vertex degrees in the 2002 snapshot of the real 
binary directed ITN (red points), and corresponding aver- 
ages over the maximum-entropy ensemble with specified out- 
degrees and in-degrees (blue solid curves), a) k]"^" versus 
fcf\ b) k™ /out versus fcf\ c) k° ut/m versus k° ut - d) k ° ut/out 
versus k° ut . 




FIG. 10: (Color online) Directed clustering coefficients ver- 
sus vertex degrees in the 2002 snapshot of the real binary di- 
rected ITN (red points), and corresponding averages over the 
maximum-entropy ensemble with specified out-degrees and 
in-degrees (blue solid curves), a) c] n versus k] n . b) c° ut ver- 
sus k° ut - c) cf c versus k™ ■ k° ut . d) cf <d versus fcj" • k° ut ■ 
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FIG. 11: (Color online) Temporal evolution of the properties 
of the total average nearest neighbor degree /j* */' * m the 
1992-2002 snapshots of the real binary directed ITN and of 
the corresponding null model with specified out-degrees and 
in-degrees. a) average of f t * ot ^ tot across all vertices (red: real 
data; blue: null model, indistinguishable from real data), b) 
standard deviation of }^° t l tot across all vertices (red: real 
data; blue: null model, overlapping with real data), c) cor- 

and k\ ot (red: real data; 
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relation coefficient between k 
blue: null model, overlapping with real data), d) correlation 



coefficient between /j' 04 / 40 ' anc i (fc* OI/tot ). The 95% confidence 
intervals of all quantities are represented as vertical bars. 
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FIG. 13: (Color online) Temporal evolution of the properties 
of the total clustering coefficient c\ ot in the 1992-2002 snap- 
shots of the real binary directed ITN and of the corresponding 
null model with specified out-degrees and in-degrees. a) av- 
erage of c* ot across all vertices (red: real data; blue: null 
model, indistinguishable from real data), b) standard de- 
viation of c\ ot across all vertices (red: real data; blue: null 
model, indistinguishable from real data), c) correlation coef- 
ficient between c* ot and k\ ot (red, upper symbols: real data; 
blue, lower symbols: null model), d) correlation coefficient 
between c| 0< and (c' ot ). The 95% confidence intervals of all 
quantities are represented as vertical bars. 
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We now consider the four directed variants & 
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k i , k i , K i .lor brevity, lor these quan- 
tities we only show the evolution of the average values, 



which are reported in Fig. 12 We find that the overall 
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behavior previously reported for the average of k i 
(Fig. [Tl^ i) is not reflected in the individual trends of 
the four directed versions of the ANND. In particular, 



^in/in 



(Fig. 12 1), k 



<J>7lf OUt 



(Fig. 12 a) and 



the averages of k i ' 

k out/out ( p. g _ 121) 

increase over a downward-shifted but 



wider range of values than that of k 



tot /tot 



whereas the 



100 

- 90 
•- SO 

70 
' 60 

50 

40 L 

1990 



1992 1994 1996 1998 2000 2002 
year 



! 120 

1 

E 115 

i uo 

S 105 



I I 1 



100 L 

1990 



1992 1994 1996 1998 2000 2002 
year 



100 | r ; 

r 90 

1 70 i 

I 60 
g 

50 

4o t , . . I : 

1990 1992 1994 1996 1998 2000 2002 
year 

d 

100 , i ' 

f 90 , 

r so 

f 70 
E 60 
50 

1990 1992 1994 1996 1998 2000 2002 
year 



FIG. 12: (Color online) Averages and their 95% confidence 
intervals (across all vertices) of the directed average nearest 
neighbor degrees in the 1992-2002 snapshots of the real binary 
directed ITN (red), and corresponding averages over the null 
model with specified out-degrees and in-degrees (blue, indis- 
tinguishable from real data), a) average offer 7 " 1 ; b) average 
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01 ' ; c) average ot k i ; d) average oi k t 
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FIG. 14: (Color online) Averages and their 95% confidence 
intervals (across all vertices) of the directed clustering coeffi- 
cients in the 1992-2002 snapshots of the real binary directed 
ITN (red), and corresponding averages over the null model 
with specified out-degrees and in-degrees (blue, indistinguish- 
able from real data), a) C*' 1 ; b) c° ut \ c) cf c ; d) c? id . 
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FIG. 15: 



(Color online) Total average nearest neighbor de- 
gree ^*°*/ tot versus total degree fc' ot in the 2002 snapshots 
of the commodity-specific (disaggregated) versions of the real 
binary directed ITN (red points), and corresponding aver- 
age over the maximum-entropy ensemble with specified out- 
degrees and in-degrees (blue solid curves), a) commodity 93; 
b) commodity 09; c) commodity 39; d) commodity 90; e) 
commodity 84; f) aggregation of the top 14 commodities (see 
Table [T] for details). From a) to f), the intensity of trade and 
level of aggregation increases. 
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FIG. 16: (Color online) Total clustering coefficient cj ot 
versus total degree fe' ot in the 2002 snapshots of the 
commodity-specific (disaggregated) versions of the real binary 
directed ITN (red points), and corresponding average over the 
maximum-entropy ensemble with specified our-degrees and 
in-degrees (blue solid curves), a) commodity 93; b) com- 
modity 09; c) commodity 39; d) commodity 90; e) commod- 
ity 84; f) aggregation of the top 14 commodities (see Table [I] 
for details). From a) to f), the intensity of trade and level of 
aggregation increases. 



average of k. 



out/in 



(Fig. 



12:) is almost constant in time. 



(Fig. 14:) display an increasing trend, whereas the 



The moderately increasing average of fc'°'/* * j s therefore 
the overall result of a combination of different trends 
followed by the underlying directed quantities, some of 
these trends being strongly increasing and some being 
almost constant. Therefore we find the important result 
that there is a substantial loss of information in passing 
from the inherently directed quantities to the undirected 
or symmetrized ones. Still, when we compare the empir- 
ical trends of the directed quantities with the random- 
ized ones, we find an almost perfect agreement. This 
implies that even the finer structure of directed correla- 
tion profiles, as well as their evolution, is reproduced in 
great detail by controlling for the local directed topolog- 
ical properties alone. 

The same analysis is shown for the total clustering co- 
and for the four directed vari- 



efficient c\ ot in Fig. 13 



ants cf\ c° ut , cf c , c™ l(i in Fig. [14} Again, we find 
that the four t emp oral trends involving the overall quan- 
tity cf* (Fig. 



131 replicate what we have found for its 
undirected counterpart c; (shown previously in Fig. [4]). 
Whe n w e consider the four inherently directed quantities 
(Fig. 14), we find that the averages of cf 1 (Fig. Il4k) and 



average ofc™ td (Fig. 141) is constant and that of c° ut 
(Fig. [Mb) is even decreasing. When aggregated, these 



different trends give rise to the constant behavior of the 
average c* ot , which is therefore not representative of the 
four underlying directed quantities. This also means that, 
similarly to what we found for the ANND, there is a sub- 
stantial loss of information in passing from the directed 
to the undirected description of the binary ITN. However, 
even the fine-level differences among the directed cluster- 
ing patterns are still completely reproduced by the null 
model. 



D. Commodity-specific binary directed networks 

We now study the binary directed ITN when disag- 
gregated (commodity-specific) representations are con- 
sidered. We repeat the analysis described above by set- 
ting A = A c with c > 0. For brevity, we report our 
analysis of the 6 commodities described in Section [ill D 
and selected from the top 14 categories listed in table [I 
(again, we found similar results for all commodities). To- 
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gether with the aggregated binary directed ITN already 
described, these commodity classes form a set of 7 differ- 
ent cases ordered by increasing trade intensity and level 
of commodity aggregation. 

In Figs . [T5| and [16] we show the behavior of the total av- 
erage nearest neighbor degree and total clustering coeffi- 
cient for the 2002 snapshots of the 6 selected commodity- 
specific networks. When compared with Figs. [7] a nd [9} 
the plots confirm what we have found in Section |III D| 
for the binary undirected case. In particular, the be- 
havior displayed by the ANND and clustering in the 
commodity-specific networks becomes less and less noisy 
as more intensely traded commodities, and higher levels 
of aggregation, are considered. Accordingly, the agree- 
ment between real and randomized networks increases, 
but the accordance is already remarkable in commodity- 
specific networks, even the sparsest and least aggregated 
ones. These results confirm that, irrespective of the level 
of commodity resolution and trade volume, the directed 
degree sequences completely characterize the topology of 
the binary directed representations of the ITN. 



V. CONCLUSIONS 

All the above results clearly imply that, in the undi- 
rected as well as the directed case, for all the years con- 
sidered, and across different commodity classes, the dis- 
assortativity and clustering profiles observed in the real 
binary ITN arise as natural outcomes rather than gen- 
uine correlations, once the local topological properties are 
fixed to their observed values. Therefore we can conclude 
that the higher-order patterns observed in all the binary 
representations of the ITN, as well as their temporal 
evolution, are completely explained by local constraints 
alone. This means that the (undirected/directed) degree 
sequence of the ITN is maximally informative, since its 
knowledge systematically conveys a full picture of the 
binary topology of the network. These results have im- 
portant consequences for economic models of trade. In 
particular, they suggest that the ITN topology should 
become one of the main focuses of international-trade 
theories. While most of the literature concerned with 
modeling international trade has focused on the problem 
of reproducing the magnitude of nonzero trade volumes 
(the most important example being gravity models [40] ) . 
much less emphasis has been put on correctly replicating 
the binary topology of the ITN, i.e. understanding the 
determinants of the process governing the creation of a 
link. However, our results clearly show that the purely 
topological structural properties (and in particular the 
degree sequence) of the ITN carry a significant amount 
of information. 

A first step in reproducing the ITN topology is the 
model in Ref. [3], where the probability pij of a trade 
relationship between two countries i and j was modeled 
as a function of the GDP values of the countries them- 
selves, and all the topological properties of the network 



were successfully replicated. Interestingly, the form of 
that function coincides with the connection probability 
of the null model considered here, shown in Appendix |B"| 
in Eq. ( B3 ) , where the role of the Lagrange multiplier Xi 



associated with fcj is played by the GDP of country i. In- 
deed, an approximately monotonic relationship between 
GDP and degree has been observed [3] , providing a con- 
nection between these two results. From another perspec- 
tive, the above remark also means that the accordance 
between the real ITN as a binary undirected network and 
its randomized counterpart is replicated under an alter- 
native null model, that controls for the empirical values 
of the GDP rather than for the degree sequence. The im- 
portance of reproducing the binary topology of trade is 
reinforced by the analysis of the ITN as a weighted net- 
work with local constraints, as we show in the following 
paper [21] . 
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Appendix A: The randomization method 

Given a real network with N vertices, there are vari- 
ous ways to generate a family of randomized variants of 
it (35] I2~4H3"3"] . The most popular one is the local rewiring 
algorithm proposed by Maslov and Sneppen [3H [32] . In 
this method, one starts with the real network and gener- 
ates a series of randomized graphs by iterating a funda- 
mental rewiring step that preserves the desired proper- 
ties. In the binary undirected case, where one wants to 
preserve the degree of every vertex, the steps are as fol- 
lows: choose two edges, say and (k,l); rewire these 
connections by swapping the end-point vertices and pro- 
ducing two new candidate edges, say and (k,j); if 
these two new edges are not already present, accept them 
and delete the initial ones. After many iterations, this 
procedure generates a randomized variant of the origi- 
nal network, and by repeating this exercise a sufficiently 
large number of times, many randomized variants are ob- 
tained. By construction, all these variants have exactly 
the same degree sequence as the real-world network, but 
are otherwise random. In the directed and/or weighted 
case, extensions of the rewiring steps defined above can 
be introduced, even if with some caution [5] [JT] . Maslov 
and Sneppen's method allows one to check whether the 
enforced properties are partially responsible for the topo- 
logical organization of the network. For instance, one can 
measure the degree correlations, or the clustering coeffi- 
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cient, across the randomized graphs and compare them 
with the empirical values measured on the real network. 
This method has been applied to various networks, in- 
cluding the Internet and protein networks [3TJ[32]. Differ- 
ent webs have been found to be affected in very different 
ways by local constraints, making the problem interesting 
and not solvable a priori. 

The main drawback of the local rewiring algorithm is 
its computational requirements. Since the method is en- 
tirely computational, and analytical expressions for its 
results are not available, one needs to explicitly gener- 
ate several randomized graphs, measure the properties 
of interest on each of them (and store their values), and 
finally perform an average. This average is an approxi- 
mation for the actual expectation value over the entire 
set of allowed graphs. In order to have a good approxima- 
tion, one needs to generate a large number M of network 
variants. Thus, the time required to analyze the impact 
of local constraints on any structural property is M times 
the time required to measure that property on the orig- 
inal network, plus the time required to perform many 
rewiring steps producing each of the M randomized net- 
works. The number of rewiring steps required to obtain a 
single randomized network is O(L) where L is the num- 
ber of links [221 EU E23, and 0(L) = O(N) for sparse 
networks while O(L) = 0(N 2 ) for dense networks [4"T] . 
Thus, if the time required to measure a given topologi- 
cal property on the original network is 0(N T ), the time 
required to measure the randomized value of the same 
property is 0(M ■ L) + 0(M ■ N T ), which is 0{M ■ N T ) 
as soon as r > 2. 

A recently proposed alternative method, which is re- 
markably faster due to its analytical character, is based 
on the maximum-likelihood estimation of maximum- 
entropy models of graphs [52]. In this method, one first 
specifies the desired set of local constraints {C a }. Sec- 
ond, one writes down the analytical expression for the 
probability P(G) that, subject to the constraints {C a }, 
maximizes the entropy 



S = -^P(G)lnP(G) 



(Al) 



G 



where G denotes a particular graph in the ensemble, 
and -P(G) is the probability of occurrence of that graph. 
This probability defines the ensemble featuring the de- 
sired properties, and being maximally random otherwise. 
Depending on the particular description adopted, the 
graphs G can be either binary or weighted, and either di- 
rected or undirected. Accordingly, the sum in Eq. (All, 



and in similar expressions shown later on, runs over all 
graphs of the type specified. The formal solution to the 
entropy maximization problem can be written in terms of 
the so-called Hamiltonian H(G), representing the energy 
(or cost) associated to a given graph G. The Hamilto- 
nian is defined as a linear combination of the specified 
constraints {C a }: 



(A2) 



where {0 a } are free parameters, acting as Lagrange multi- 
pliers controlling the expected values {(C a )} of the con- 
straints across the ensemble. The notation C Q (G) de- 
notes the particular value of the quantity C a when the 
latter is measured on the graph G. In terms of H(G), the 
maximum-entropy graph probability P{G) can be shown 
to be 



P(G) 



-H(G) 



(A3) 



where the normalizing quantity Z is the partition func- 
tion, defined as 



-H(G) 



(A4) 



Third, one maximizes the likelihood P(G*) to obtain the 
particular graph G* , which is the real- world network that 
one wants to randomize. This steps fixes the values of 
the Lagrange multipliers that finally allow to obtain the 
numerical values of the expected topological properties 
averaged over the randomized ensemble of graphs. The 
particular values of the parameters {0 a } that enforce the 
local constraints, as observed on the particular real net- 
work G*, are found by maximizing the log-likelihood 



A = lnP(G* 



-H{G*)-\nZ 



(A5) 



to obtain the real network G*. It can be shown [TH] that 
this is equivalent to the requirement that the ensemble 
average (C a ) of each constraint C a equals the empirical 
value measured on the real network: 



(C a )=C a (G*) Va 



(A6) 



Note that we generally adopted a simplified notation by 
writing C* (or even only C a ) instead of C a (G*) for the 
empirically observed values of the constraints (see for in- 
stance Sections II A and II B I. Once the parameter values 



are found, they are inserted into the formal expressions 
yielding the expected value 



(X) = ^X(G)P(G) 



(A7) 



G 



of any (higher-order) property of interest X. The quan- 
tity (X) represents the average value of the property X 
across the ensemble of random graphs with the same av- 
erage (across the ensemble itself) constraints as the real 
network. In what follows we provide a detailed account of 
the expressions for the randomized properties appearing 
in our analysis. 

Technically, while the local rewiring algorithm gen- 
erates a microcanonical ensemble of graphs, contain- 
ing only those graphs for which the value of each con- 
straint C a is exactly equal to the observed value C a (G*), 
the maximum-likelihood method generates an expanded 
grandcanonical ensemble where all possible graphs with 
N vertices are present, but where the ensemble aver- 
age of each constraint C a is equal to the observed value 
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C a (G*). One can show that the microcanonical approach 
converges to the grandcanonical one as the number of 
microcanonical randomization steps increases [22 . How- 
ever, the maximum-likelihood one is significantly faster. 
Importantly, enforcing only local constraints implies that 
-P(G) factorizcs as a simple product over pairs of vertices. 
This has the nice consequence that the expression for (X) 
is generally only as complicated as that for X. In other 
words, after the preliminary maximum-likelihood estima- 
tion of the parameters {0 a }, in this method the time re- 
quired to obtain the exact expectation value of an 0(N T ) 
property across the entire randomized graph ensemble is 
the same as that required to measure the same property 
on the original real network, i.e. still 0(N T ). Therefore, 
as compared to the local rewiring algorithm, which re- 
quires a time 0(M-N T ), the maximum-likelihood method 
is 0(M) times faster, for arbitrarily large M. 



Appendix B: Binary undirected properties 

In the binary undirected case, each graph G is com- 
pletely specified by its (symmetric) Boolean adjacency 
matrix A. The randomization method described above 
proceeds by 

1. specifying the degree sequence as the constraint: 
{Co] — {ki}. The Hamiltonian therefore reads 



H(A) = *ifci(A) = ]T + 9i)t 



(Bl) 



i j<i 



and one can show (32] that this allows to write the 
graph probability as 



P(A) mi/. ;l „,,:' 



where 



Pij = 



X i X j 

1 I tJC 2 Xj -? 



(B2) 



(B3) 



(with Xi = e~ 0i ) is the probability that a link 
exists between vertices i and j in the maximum- 
entropy ensemble of binary undirected graphs, sub- 
ject to specifying a given degree sequence as the 
constraint; 

2. solving the maximum-likelihood equations, by set- 
ting the parameters {x{\ to the values that maxi- 
mize the likelihood P(A*) to obtain the real net- 
work A* [TBI [22] . These values can be found as the 
solution to the following set of N coupled nonlinear 
equations [T8] : 



J^=fc. ( (A*) 

I X i X j 



Vi 



(B4) 



where {fcj(A*)} is the empirical degree sequence of 
the real network A*. For a detailed analysis about 



solving such system see (for a discussion about 
the existence of solutions) and [22] (for a discussion 
about the convergence of the algorithm) . In princi- 
ple, dimensionality and memory problems can arise 
when N is too large (luckily, this is not the case of 
the ITN considered here). In such a case, the sys- 
tem can be re-written with a lower number of equa- 
tion to solve. In fact, the hidden variables of the 
vertices with the same degree have the same value. 
So, one can straightforwardly solve the system only 
for them [Tgj ; 

computing the probability coefficients , by insert- 
ing the values {xi} into eq. (B3) which allows to 



easily compute the expectation value (X) of any 
topological property X analytically, without gen- 
erating the randomized networks explicitly |22j . 
With this choice, eq. (B3) yields the exact value 



of the connection probability in the ensemble of 
randomized networks with the same average de- 
gree sequence as the empirical one. Note that p^ is 
the probability of a link between vertex i and ver- 
tex j in the grandcanonical ensemble (which is di- 
rectly obtained analytically) , and not the frequency 
of such a link in the corresponding microcanonical 
ensemble (which would require the explicit genera- 
tion of artificially rewired networks) . In ref . [22j , it 
was shown that the microcanonical frequency con- 
verges to pij asymptotically as the number of ran- 
domization steps in the microcanonical algorithm 
increases. Equation (B4| shows that, by construc- 
tion, the degrees of all vertices are special local 
quantities whose expected and empirical values are 
exactly equal: (ki) — ki. It follows that the p^ 
coefficients can be calculated by using any of the 
networks in the corresponding microcanonical en- 
semble with constrained degree sequence: the ex- 
pected values of the high-order properties will be 
the same; 

4. computing the expectation values of higher-order 
topological properties, as in Tabic [Tl] The expres- 
sions are derived exploiting the fact that (a^) = 
Pij, and that different pairs of vertices are statisti- 
cally independent, which implies (aijOki) — PijPkl 
if (i — j) and (fc — I) are distinct pairs of vertices, 
whereas (a i:j a kl ) = (a?-) = (a tJ ) = p^ if (i - j) 
and (k — I) are the same pair of vertices. Also, the 
expected value of the ratio of two quantities is ap- 
proximated with the ratio of the expected values: 
(n/d) » (n)/(d). 



Appendix C: Binary directed properties 

In the binary directed case, the above results can be 
generalized as follows. Each graph G is completely spec- 
ified by its Boolean adjacency matrix A, which now is 
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Empirical undirected properties 


Expected undirected properties 
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Empirical directed properties 


Expected directed properties 
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j out/out X^j'^i a ijkj Ut 
i - kf* 


n out/ out \ Y J j 1 LiVijk° ut 
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j tot/tot Y.j^i( a ij+ a ji) kt j 0t 
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^ 1 ' Sfc^i.j PjiPki 
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mid l^j^i l^k^i.j a ik a jk a ji 


1 mid\ 'iZj^LiY.kjti.jPikPjkPji 


C l - fc in fc out_ fe « 




c tot _ 2_'j^i Ufc^iJ ( a ij+ a jiH a jfc+ a fcjH a fci+ a ifc) 




2 feJ ot (fcJ ot -l)-2fc1^ 


2 EjVi 'Ek^i,j(PjiPki+PijPik) + 2( k l nk i Ut )-2'Ej^i PijPji 



TABLE II: Expressions for the empirical and expected properties in the binary (undirected and directed) representations of 
the network. 



in general not symmetric. The maximum-likelihood ran- 
domization method [22 proceeds in this case by 

1. specifying both the in-degree and the out-degree se- 
quences as the constraints: {C a } = {fc™,fc° n *}. 
The Hamiltonian takes the form 



H(A) 



E 



[0f fcf (A) + 8° ut k° ut (A)] 



- ££(« 



+ Vj )ai 



(Cl) 



The above choice leads to the graph probability [22] 



/'A, I J 1 J /> 1 



1 — Oi 



where 



CtVj 



1 



(C2) 



(C3) 



(with xi = e~ e »°" and yt = e" 9 *") is the probabil- 
ity that a link exists from vertex i to vertex j in 
the maximum-entropy ensemble of binary directed 
graphs with specified in- and out-degree sequences; 

2. solving the maximum-likelihood equations, by set- 
ting the parameters {x^ and {yi} to the values 



that maximize the likelihood P(A*) to obtain the 
real network A* [18 , 22J . These values can be found 
as the solution to the following set of 2N coupled 
nonlinear equations [15] : 



£ , •''' // \. *r';A / (CM) 



(fef) ^ Y : — = fcj w (A*) Vt (C5) 

^ 1 + Sj-W 



where {fc-"(A*)} and {fcf"'(A*)} are the empirical 
degree sequences of the real network A*. Again, 
for a detailed analysis about solving such system 
see [HES]; 

computing the probability coefficients pij, by insert- 
ing the values {x^ and {yi} into Eq. (C3) which 



allows to easily compute the expectation value (X) 
of any topological property X analytically, with- 
out generating the randomized networks explicitly 
[2"2"] . So Eq. (C3) yields the exact value of the con- 



nection probability in the ensemble of randomized 
directed graphs with the same average degree se- 
quences as the empirical ones and Eqs. (C4)-(C5) 



show that, by construction, the in-degrees and out- 
degrees of all vertices are special local quantities 



18 



whose expected and empirical values are exactly 
equal: (fe| n ) = fcf" and (k° ut ) = k° ut . It follows 
that the pij coefficients can be calculated by using 
any of the networks in the corresponding micro- 
canonical ensemble with constrained in-degree and 
out-degree sequences: the expected values of the 
high-order properties will be the same; 



4. computing the expectation values of the higher- 
order topological properties, as in Table [TTJ by using 
the same prescription as in the undirected case plus 
the additional care that now (i — j) and (j — i) are 
different (and statistically independent) directed 
pairs of vertices. Therefore (ayajj) — PijPji- 
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